function tstat = glmtstat(c, b, iXX, se)
% tstat  - compute GLM t statistic
%
% FORMAT:       tstat = glmtstat(c, b, iXX, se)
%
% Input fields:
%
%       c       1xP contrast vector (right padded with zeros)
%       b       beta estimates (double or single)
%       iXX     pinv of X'*X
%       se      standard error of estimation (double or single)
%
% Output fields:
%
%       tstat   t-statistic
%
% See also @BVQXfile/private/sdm_CalcBetas

% Version:  v0.7f
% Build:    8110521
% Date:     Nov-05 2008, 9:00 PM CET
% Author:   Jochen Weber, SCAN Unit, Columbia University, NYC, NY, USA
% URL/Info: http://wiki.brainvoyager.com/BVQXtools

% argument check
if nargin ~= 4 || ...
   ~isa(c, 'double') || ...
   (~isa(b, 'double') && ...
    ~isa(b, 'single')) || ...
   ~isa(iXX, 'double') || ...
   (~isa(se, 'double') && ...
    ~isa(se, 'single')) || ...
    isempty(c) || ...
    isempty(b) || ...
    ndims(iXX) ~= 2 || ...
    size(iXX, 1) ~= size(iXX, 2) || ...
   ~any(size(b) == size(iXX, 1)) || ...
    numel(se) ~= (numel(b) / size(iXX, 1))
    error( ...
        'BVQXtools:BadArgument', ...
        'Missing, invalid or bad sized argument given.' ...
    );
end

% check sizes
csize = size(iXX, 1);
c = c(:)';
if numel(c) < csize
    c = [c, zeros(1, csize - numel(c))];
else
    c = c(1:csize);
end

% determine beta dimension
bd = find(size(b) == csize);

% if needed, reorder beta
if bd(1) ~= 1
    b = reshape(permute(b, [bd(1), setdiff(1:numel(size(b)), bd(1))]), ...
        [csize, numel(se)]);
end

% calc stats
tstat = reshape((c * b) ./ (se(:)' * sqrt(c * iXX * c')), size(se));
